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ABSTRACT. Global convective flows in the solar convection zone have been 
predicted by theoretical interpretations of the global-scale ordering of 
magnetic fields and activity centers and by theoretical analyses of 
rotating convection zones. Direct evidence of these flows in the photo- 
sphere has not previously been found despite several long-term efforts. 
The signatures of such flows have now been detected by analyzing the 
daily series of low-resolution Dopplergrams obtained at the Wilcox Solar 
Observatory at Stanford University. The signatures are patterns of 
alternating east and west flows with amplitudes on the order of 25 m/s 
and longitudinal extent of about 30 degrees. The patterns move across 
the disc at approximately the solar rotation rate and have lifetimes of 
at least several rotations. Boundaries of the fast and slow flows are 
often associated with large magnetic active regions. 


The detection and study of global convective wave flows on the sun is of 
great importance for understanding the basic mechanisms of the genera- 
tion of magnetic fields and magnetic activity on the sun (Yoshimura, 
1971). They have been considered as the most likely kind of fluid 
motions to drive the solar dynamo which generates magnetic fields on the 
sun (Yoshimura, 1972, Yoshimura, 1983). This global convection has 
escaped detection despite a number of long-term observing and analysis 
programs (Durney et al., 1985, Howard and Yoshimura, 1976, Howard and 
LaBonte, 1980, LaBonte et al., 1981, Snodgrass and Howard, 1984). 

Observations made at the Wilcox Solar Observatory at Stanford 
University have been examined for the signature of global scale convec- 
tive motions. This series of observations began in 1976 and has been 
found to have lower noise in the detection of solar velocities than the 
instruments used in previous analyses (Scherrer et al., 1980). We now 
have a sufficient time span of observations to study the evidence for 
global-scale horizontal motions superimposed on solar rotation. The 
observations are Doppler-shift measurements which yield low-resolution 
full-disc daily velocity maps. Clear velocity patterns with amplitudes 
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of at least 25 m/s have been found in zones near the equator. These 
patterns can be seen to move across the disc with solar rotation. They 
have a large latitudinal extent, crossing the equator, but are not 
aligned on meridians. There are from 4 to 7 alternating patterns around 
the sun (Scherrer and Yoshimura, 1985, Yoshimura et al., 3-15). 

There are several fundamental difficulties in detecting global- 
scale small-amplitude horizontal motions on the sun. One problem is 
caused by the large contributions of solar rotation (about 2000 m/s at 
the equator) and the limb shift (a "W'-shaped red-shift profile across 
the disc with an amplitude of 600 m/s at the limb and 30 m/s at the disk 
center) . Other sources of solar "noise" are the motions of 5-minute 
oscillations and supergranulation. These noise velocities are larger in 
amplitude but smaller in scale than the postulated horizontal flows. 
The large aperture used for the Stanford observations greatly attenuates 
these signals. A third solar source of velocity noise comes from the 
red-shifts associated with regions of intense magnetic fields. These 
"downflows" have been examined by a number of authors (Cavallini et al., 
1985, Labonte and Howard, 1982) and are generally believed to be due to 
a local reduction in small scale convective motions in the vicinity of 
magnetic fields. 

A final difficulty is that we are restricted to measuring the 
line-of-sight component of the motion. In order to infer individual 
velocity components, one must assume either spatial continuity or tem- 
poral continuity during the disc passage of individual features. In the 
present work, we have used two different methods in an attempt to under- 
stand and avoid the weaknesses of each. In both approaches, we assume 
the rotation and limb shift signals are constant over several rotations, 
and that the remaining horizontal motions are organized on a scale much 
larger than supergranulation. In the first approach (Scherrer and 
Yoshimura, 1985) we further restrict ourselves to horizontal velocity 
structures with a large east-west extent. We assume that by averaging 
over a large window (up to a solar radius) the components other than 
east-west motions will cancel. In the second approach (Yoshimura et 
al., 3-15), we do not assume such a large spatial extent, but we assume 
that the structures are stable during their disc passage. The two 
methods have different sensitivities to noise and changing solar velo- 
city fields. We refer to the two methods as the slope method and the 
decomposition method. 

The Stanford velocity maps are a product of a 9-year synoptic 
series of observations. The observations and instrument have been 
described previously (Duvall, 1979, Scherrer et al., 1980) and will be 
only briefly described here. The instrument is a scanning magnetograph 
attached to a high dispersion spectrograph. The raw data consist of 
intensity, magnetic field, spectral line position, position on the solar 
image, and time. Daily scans are made with a - three-arc-minute square 
aperture. The solar image is scanned in 21 steps of 90 arc-sec in the 
solar east-west direction and stepped to 11 scan lines in the north- 
south direction. In this way we obtain simultaneous magnetograms and 
Dopplergrams on every suitable day. The observations are made in the 
absorption line of neutral iron near 525.02nm. The calibration of the 
line shift signal in m/s and removal of the earth's motion is done as 
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described previously. After excluding observations contaminated by 
clouds we have 1494 daily Dopplergrams for the present study. 

Previous analyses have removed the large rotation and limb-shift 
signal by fitting and subtracting a particular functional form of dif- 
ferential rotation and limb-shift from each Dopplergram. To avoid any 
spurious signal from an improper fit and to avoid the possibility of 
removing some of the desired signal as rotation, we have adopted a new 
procedure. The first step is to subtract the 54-day running mean of the 
observed velocity at each point on the disc. This procedure removes any 
signal dependent on disc position, such as rotation and limb shift, 
without making any assumptions regarding its functional form. Next, we 
subtract from each scan point residual the mean of all residuals within 
the central 80$ of that Dopplergram. This step is necessary because we 
have no absolute velocity reference for our measurements. To reduce the 
contamination of our horizontal velocity signal by active-region red- 
shifts, we next discard all scan points which have a velocity residual 
amplitude exceeding 30 m/s or a magnetic field exceeding 20 gauss. An 
average of 27 percent of the scan points were discarded from each Dop- 
plergram. Discarding these points may remove evidence of the larger 
amplitude flows, but we believe it is necessary to remove most of the 
contamination from the active regions. 

The resulting series of residual Dopplergrams is then analyzed for 
horizontal velocity components by fitting the east-west component of 
motion. In the first approach (slope method) we fit the east-west com- 
ponent in a restricted portion of the solar disc with a large longitudi- 
nal window at a fixed latitude. This is done by finding (by least 
squares) the east-west component of velocity (V west ) such that: 
v obs = V west sin(L)cos(B 0 ) where L is the longitude measured from disc 
center and B 0 is the latitude of disc center. The regions included in 
each fit were east-west strips 810 arc-seconds long by 180 arc-seconds 
wide. These fits were computed for the nine central scan lines and all 
possible east-west positions for each Dopplergram. From these fits we 
can display daily maps of east-west motions. The relatively large area 
used in these fits greatly attenuates the noise from oscillations and 
supergranulation. It also acts as a spatial filter allowing only very 
large-scale structures to produce a signal. The method is not sensitive 
to day-to-day variability in instrumental zero or noise from changing 
solar velocities, but it is sensitive to contamination from line shifts 
not associated with east-west motion. 

We also derived east-west motions by a second approach (decomposi- 
tion method) in which we examine the Doppler-shift variation of each 
location on the sun as it crosses the disc. This is accomplished by 
using all Dopplergrams having a given Carrington coordinate within 0.8 
R 0 of disc center to find the inferred components 1), V, and W that best 
fit (least squares) the observed velocity residual V 0 b s : 
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Figure 1: Velocity residual maps for September 19 to 23, 1979. The top 

row is the line-of-sight "velocity" with limb-shift and rotation removed 

(positive away from observer). The middle row if the computed east-west 

velocity (positive to the west) and the bottom row is the magnetic 

field. The contour levels are 10 m/s for the velocity maps and 0, +.1, • 

2, 5, 10 ... gaus for the matnetograms . The regions excluded from the 

analysis are shaded on the magnetograms. 


v obs = 0cos(B o )sin(L) 


+V(cos(B 0 )sin(B)cos(L) - sin(B 0 )sin(B)) 

-W(cos(B 0 )cos(B)cos(L) + sin(B 0 )cos(B)) 

where B is the latitude and U, V, W are the westward, northward, and 
upward components respectively. This method has the advantage that it 
cleanly separates the westward velocity from active-region "downdrafts" 
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but has the disadvantage that it assumes a static velocity structure, no 
day-to-day zero level errors, and the Carrington rotation rate. 

Using the slope method, Scherrer and Yoshimura (1985) have shown 
that there are global-scale velocity features that are stable for their 
transit across the disk. This is the first requirement for detecting 
large scale flows. Figure 1 shows data from 5 consecutive days, Sep- 
tember 19 to 23, 1981. The top row shows residual line-of-sight signals 
(before excluding points with large residuals or fields), the center row 
shows V wes t as described above, and the bottom row shows the correspond- 
ing magnetograms for comparison. The V wes t signal can easily be seen to 
move across the disc as the sun rotates. The V wes t signal is shown as 
positive values with solid contours in regions of westward motion (i.e. 
motion in the same direction as rotation) . 

To study the morphology of the horizontal flows we have con- 
structed synoptic charts of both V wes t and U in the same way as is usu- 
ally done for magnetic fields. These charts show a map of the solar 
surface made by averaging the signals from each day's observations for 
each latitude-longitude bin. The synoptic charts are computed on a grid 
of 5-degrees in longitude by 30 equal steps in sine-latitude. Some 
additional temporal and spatial averaging results from this process. To 
reduce the small-scale noise a further averaging with a 40-degree-in- 
longitude low-pass filter was applied to each latitude strip. 

Figure 2a shows V we3t for three consecutive rotations from the 
fall of 1979, Figure 2b shows three rotations from the fall of 1981, and 
Figure 2c shows the fall of 1984. The velocity pattern shows 4 or 5 
structures of alternating eastward and westward motion per rotation. 
These structures are most clearly seen in low northern latitudes. They 
can be followed for several rotations. 

Synoptic maps computed with the decomposition method have roughly 
the same character but differ in detail. They differ most in the vicin- 
ity of active regions, where the slope method could be contaminated by 
"downflows" if their spatial extend is large, and the decomposition 
method would be contaminated if they are changing. 

This signal appears to persist for more than one rotation. Flows 
of about the same magnitude are generally present during the past six 
years during periods of both intense and low magnetic activity. For 
instance, the signal is present in the fall of 1984 when activity was 
quite low indicating that the while there is a correlation with the mag- 
netic structures, we really are measuring different physical quantities. 
A full understanding of the relation between the flow patterns and mag- 
netic activity will require more analysis and more complete observations 


Figure 2: Synoptic maps of east-west velocity V we3 t for three consecu- 
tive rotations in 1979 (part a), in 1981 (part b), and in 1984 (part c). 
Westward moving regions are shown in unshaded. Eastward moving regions 
are shaded dark. The contour levels are ±10, 20, 50, ... m/s. 
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sequences than are presently available. 

The intervals presented in Figure 2 were chosen both because of 
sufficient data coverage to show nearly complete maps and because an 
initial survey had identified these intervals as times of simple modal 
structure. This survey was done by computing power spectra of V wes t for 
intervals of nearly continuous data coverage. It is interesting to note 
that the day-to-day variation in the whole-disc fit equatorial rotation 
rate (Scherrer et al., 1980) is very similar to the meridional average 
of V west for the intervals examined. This was examined both by compar- 
ing power spectra and by inspection of the raw data. Essentially all of 
the day-to-day variation in the full disc fit of rotation comes from the 
velocity structures described here. 

These velocity patterns have a large amplitude compared to the 
upper limits reported in previous analyses of observations from Mt. Wil- 
son (Snodgrass and Howard, 1984). However, we note that the present 
analysis of Stanford observations finds horizontal flow patterns which 
are different in structure, and larger in scale than those searched for 
in the Mt. Wilson studies. 

The complete understanding of the flow patterns reported here will 
require an extended set of high quality continuous observations 
throughout a solar cycle. The present results have been obtained with 
the instrument at Stanford. This instrument at present is limited to 
low resolution measurements and seldom has more than 30 or 40 days 
without interruption of a day or two. We are taking steps to increase 
both the resolution and the Doppler zero stability for future observa- 
tions. We are also examining the comparable data of improved quality 
available from Mt. Wilson since 1982. We hope to be able to test our 
interpretation of the Stanford observations as signatures of the global 
convection and to begin the process of examining the morphology of these 
flows by comparisons with other data and further observations. 
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Figure 1 

The power spectrum of the full data, from 1000 pHz to 5000 pHz 


5000 



Figure 2 


The window function: power spectrum of a 3000 jjHz sine wave put into 
the actual data window. 


to the main peaks greatly affect the fine scale structure of the spec- 
trum and make the identification of individual peaks difficult. 

In addition to the complications introduced by the window, the 
effects of splitting due to rotation and lifetimes, each producing mul- 
tiple peaks with sidelobes, as well as noise make the spectrum very com- 
plex. Several computerized deconvolution and peak-subtracting tech- 
niques failed to make significant inroads on the problem (Scherrer 
1985). A scheme combining computed power spectra and manual pattern 
recognition was decided upon as the best method to identify real solar 
peaks. Figures 3 and 4 are typical examples of the graphs used. These 
figures each show a section of the spectrum computed with several reso- 
lutions. Frequency runs along the abscissa, with each panel covering 
the same span, from 1432 jjHz to 1572 pHz in figure 3, and from 3200 pHz 
to 3340 pHz in figure 4. We realized that the full resolution spectrum 
(panel E) is too complex to enable recognition of any given peak as a 
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Figure 3 

Plot of data in the frequency range from 1432 jjHz to 1572 jiHz : 

A) average power spectrum of 4 day chunks of data, B) peak search on A, 
C) average power spectrum of 8 day chunks of data, D) peak search on C, 
E) full resolution power spectrum, F) peak search on E. 


solar oscillation mode. At lower resolution, however, it becomes much 
clearer to see where there is power in the spectrum. Thus, the average 
of 6 spectra of 4-day segments of the data is plotted in panel A. Panel 
C shows the average of 3 spectra of 8-day segments of the data. To 
further aid in locating modes a simple computer program was used to scan 
the spectra, looking for structures which correspond to the main 
features of the window function, i.e. a central peak with two sidelobes 
at 1/day distance and of about 50$ magnitude. The scan allowed the 
sidelobes to be as large as 90$ or as small as 20$ of the central peak. 
The program took the value of the spectrum at a frequency if it 
corresponded to a possible central peak and zero otherwise. The results 
of this search are plotted below the corresponding spectrum, in panels 
B, D and F. This composite graph proves to be a very useful tool in 
determining the frequencies of the modes. The first four panels, A-D, 
serve to locate and emphasize small frequency ranges where there is 
power. Using this information and a transparency with the window func- 
tion plotted on it, one can then pinpoint the group of peaks (and its 
sidelobes) in the high- resolution full spectrum, panels E and F, that 
correspond to a solar p-mode. The frequency recorded for each mode was 



Figure 4 

As figure 3, but for the range from 3200 jjHz to 3340 pHz. 


that of the central peak or approximate center of the group of peaks. A 
conservative estimate of the error in the frequency of any given mode is 
2 pHz. 

Analyzing the spectrum, we located modes ranging from 1035 pHz up 
to 5008 pHz. Identification of the modes was accomplished using previ- 
ously published frequencies (Duvall and Harvey 1983, Duvall - private 
communication) for modes in the range from about 2000 to 4000 pHz as 
"anchors" for the determination of order n and degree 1. The fact that 
the Stanford instrument is most sensitive to modes of degree 2 through 5 
prevents confusion with higher degrees, although those modes certainly 
contribute to the "noise". Modes with frequencies outside of this range 
were then identified by making an echelle diagram (figure 5) and follow- 
ing the curve described by different n values of a given 1. It is 
important to note that although there certainly is some bias in the 
location and identification of the modes with previously published fre- 
quencies, there is no such bias in the frequencies of modes of lower and 
higher order, since there was no preconception of where they should be. 
Table I lists the frequencies of the identified modes. 

In addition to identification of the modes, we have also examined 
the data for evidence of rotational splitting. As mentioned before, the 
power of each mode is split into a multitude of peaks resulting from a 
combination of mode lifetimes, window function, and noise, as well as 
rotational splitting. This confusion effectively prevents one from 
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Figure 5 

Echelle diagram of the identified modes, each strip covers 136 pHz, 
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1691 

10 
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1726 
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11 

1802 

1859 

1916 


12 


2001 

2052 

2105 

13 


2139 

2187 

2236 

14 

2217 

2274 

2052 

2372 

15 

2352 

2411 

2457 

2505 

16 

2486 

2543 

2591 

2640 

17 

2621 

2677 

2727 


18 

2754 

2812 

2864 


19 

2889 

2948 

3001 

EmSi 


Table I - Frequencies of identified modes 


assigning m-values to given peaks and directly measuring the splitting. 
Rather, we have analyzed the spectrum using autocorrelations. Using the 
determined frequencies we calculated the autocorrelation of a small sec- 
tion of the full spectrum around each mode. All autocorrelations with a 
given degree 1, but different order n, were then averaged together. The 
assumption is that any rotational splitting will not change much with 
order n. The results are plotted in the first four panels of figure 6. 
An autocorrelation of the window function is shown in the last panel of 
figure 6. It can be seen that the window is responsible for much of the 
systematic splitting, as evidenced by the peaks in its autocorrelation 
which are repeated and actually dominate the autocorrelations of the 
modes. However, all degrees show a peak around 0.8 pHz, which is close 
to twice the rotational splitting as measured by Duvall and Harvey 
(1984). Since the Stanford instrument is sensitive only to modes with 
1+m even, this peak can be identified as the result of rotational split- 
ting. 


2. CONCLUSION 

We have determined the frequencies of over 20 modes which have not been 
identified before. The majority of these modes lies in the frequency 
range from 1035 pHz to 1800 pHz. The echelle diagram of the modes shows 
clearly the change in slope at low n, as the modes become spread further 
apart than 136 pHz. Also, we have been able to see a characteristic 
splitting which is consistent with the rotational splittings measured by 
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Figure 6 

Autocorrelations of the power spectrum around different degrees 1 and 
the window function. 


Duvall and Harvey (1984). Future examination of other Stanford data 
will hopefully refine the values, although no other year has had data as 
complete as 1984. 
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ABSTRACT. Current observational questions in asteroseismology require 
high spectral resolution that can only be obtained with observations 
spanning many days or months. The primary constraint in the full utili- 
zation of single mid-latitude observing sites is the presence of diurnal 
data gaps. Several methods for removing the effect of these gaps in the 
spectra obtained from velocity observations have been suggested. In the 
case of data coverage of less than 50 % , none of the methods considered 
has been successful at unambiguously recovering the true spectrum from 
test spectra with realistic complexity. The limitations of these 
methods and their applicability to the helioseismology problem is dis- 
cussed. 


The presence of data gaps is particularly troublesome when computing 
power spectra of multi-day observations of solar or stellar oscilla- 
tions. If the expected lifetime of the observed modes is longer than 
individual uninterrupted observing sequences, one naturally wants to 
combine many days of data to obtain higher spectral resolution. The 
gaps from nights or poor weather produce a complicated structure of 
sidebands in the power spectrum unless the gaps are filled with 
estimated data. Since the power spectrum is computed from the Fourier 
transform, and the observed transform can be considered as the convolu- 
tion of the true spectrum and the transform of the observing window, 
filling the gaps is essentially a deconvolution process. Since there 
are an unlimited number of ways to fill the gaps with data, most of them 
non-physical, we must have some guidance from our a priori knowledge of 
the sun or star in question and we must have some method for incorporat- 
ing this knowledge in an analytical procedure. Thus the whole problem 
of deconvolution or gap-filling is one of explicitly recognizing the 
assumptions about the star and finding objective ways to incorporate 
those assumptions in the analysis. 

This problem is not unique to asteroseismology. In fact, there is 
a significant body of literature dealing with the deconvolution problem. 
The specific problems of asteroseismology are different from many of the 
previous applications however, so care must be exercised in applying 
previously developed methods. Some of the particular problems of solar 
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and stellar observations are low filling factors (i.e. large gaps rela- 
tive to intervals of continuous data) , very complex spectra with 
thousands of modes allowed over several decades in frequency, line 
separation on the same order as the main contributions from typical 
observing windows , and low signal to noise in frequency bands that con- 
tain important oscillation modes. A number of gap-filling methods have 
been suggested for the asteroseismology problem. These include maximum 
entropy methods of "predicting" the data in the gaps (Fahlman and 
Ulrych, 1982), variations of the CLEAN algorithm of iterative peak sub- 
traction (Schemer, 1984), rearranging the data to fill the gaps at the 
expense of resolution (Kuhn, 1982), and iterative multiplicative decon- 
volution (Connes and Connes, 1984). 

Brown and Christensen-Dalsgaard (1985) have completed a thorough 
study of the applicability of the maximum entropy method to the problem. 
They found that complex spectra can be correctly recovered provided the 
filling factor is high (greater than 0.8) and the signal-to-noise ratio 
is high (greater than 100). Fortunately, these conditions are expected 
to be met for p-mode observations from an earth-bound network They found 
that as the filling factor dropped below 0.5 or the s/n dropped below 
about 30 the error rate increased rapidly. Kuhn found similar results 
when applying the data-swapping method to data with a low filling factor 
(Kuhn et al. , 1985) . 

Delache and Schemer (1983) successfully applied the CLEAN-like 
iterative peak subtraction method to low-filling factor observations in 
the long-period spectral region but that method was not successful when 
applied to p-mode spectra (Henning and Schemer, 1985). Experiments 
with linear iterative deconvolution methods have not been as successful 
as the CLEAN technique. There is a tendency for noise peaks to grow 
rather than to be simply extracted unchanged as in the CLEAN method. 
This is a characteristic of linear iterative methods (Jansson, 1984). 

In order to push these procedures to the regime of low-filling 
factor, the methodology of inserting a. priori knowledge into the 
analysis must be explicitly examined. There are different assumptions 
implicitly made by the various methods. The maximum entropy method 
assumes that the signal is locally stationary. This is a fairly weak 
assumption and leads to failure when the gaps are large. The stronger 
assumption of completely stationary signals made in the CLEAN method can 
do a better job of recovering the true spectrum with low filling fac- 
tors, but only in cases where the assumption is valid. Thus, the long 
period success and the p-mode failure of the method. 

Particularly strong claims have been made for the advantages of 
the non-linear spectral deconvolution method (Connes and Connes, 1984). 
This has lead to a careful examination of the applicability of that 
method to the problem. The method can be summarized as the iterative 
task of finding a spectrum that when convolved with the window function 
yields the observed spectrum. The iteration is done under some con- 
straints that provide the mechanism for including prior knowledge. This 
method has the advantage of explicitly specifying the assumptions. It 
has been applied for some time to the deconvolution of instrumental pro- 
files from spectra and has been described by (Blass and Halsey, 1981) 
for that application. The general method has been examined by a number 


of authors in other applications and it has been shown that convergence 
can be guaranteed in many cases (Sanz and Huang, 1983, Schafer et al., 
1981). A general discussions of both subtractive and multiplicative 
deconvolution techniques have been presented by Bates and colleagues 
(Bates et al., 1982a, 1982b, 1982c) and by Jansson (1984). 

Connes and Connes suggested applying the multiplicative deconvolu- 
tion method independently to the real and imaginary parts of the Fourier 
transform of oscillation observations. It turns out that constraints 
suggested do not quite match the problem. In particular, the method 
assumes that the true signal is band-limited, and that power outside the 
known band comes from the gaps. Unfortunately, it is not possible to 
arrange the observations of the sun (or many stars) to match this cri- 
terion. The p-mode frequency band, for example, extends from possibly 
300 micro-Hz through 4500 micro-Hz. The sidebands from nightly gaps are 
spaced at about 11.6 micro-Hz. Thus most of the sidebands are within 
the desired frequency band. The method probably works but is expected 
to be rather inefficient. Thus, this method at present does not solve 
our problem (J. Connes, private communication, 1985). 

What has yet to be done, is to find a methodology for tuning an 
analysis to specific prior knowledge or assumptions. For example, if we 
were asking what is the rotational splitting for solar p-modes, then 
since we know the central frequencies for each multiplet (to a few 
micro-Hz) we should be able to tell our analysis procedure just where 
power is from gaps and where it is probably of solar origin. If we 
believe that the lifetime of low-frequency p-modes is several months but 
the lifetime of high-frequency modes may be only a few days, we should 
be able to find peak centroids using varying resolution across the p- 
mode band. Perhaps the iterative methods can be adapted to such more 
problem-specific constraints. In any case, we must always be aware of 
the explicit and implicit assumptions in our procedures. In summary, 
two things are clear: First, if we do not attempt to account for the 
sidebands introduced by the data gaps we are making the implicit assump- 
tion that the sun does not oscillate at night and second, we do not yet 
have a reliable, objective method for filling the gaps in the case of 
low filling factor. 
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